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Abstract 

The coupled multi-physics model for offshore wind turbine 
system has been developed and presented in this paper. The 
subsystem models are developed initially separately and 
coupled together based on the partitioned method. The 
multi-physics model consists of four subsystems, structural 
dynamics of turbine, described by finite element method, 
aero dynamics of wind which is modelled based on blade 
element momentum and dynamics stall model, hydro- 
dynamics of wave which is modelled based on boundary 
element method and soil dynamics of sea bed is modelled 
based on a coupled of finite element and scaled boundary 
finite element. The software components communicate with 
the so called Component Template Library (CTL). The 
simulation of entire model of offshore wind turbine is 
performed in time domain and results are presented at the 
end of this paper. 
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Introduction 

Due to the increasing requirement for energy over the 
world, humans always seek new sources of energy 
having no impact to environment. The wind energy 
technology is a renewable energy technology which 
has less impact on environment. In history, wind 
turbine was used only for agriculture. As people are 
concerning and searching for new sources of energy, 
the development of wind turbine for running the 
electrical power generator is initialized. In the long 
period of research and development, wind turbine is 
relatively small in the past compared to present 
turbine. Previously, there were turbines with rotor 
diameter of 10 meters; while presently they are very 
large and more flexible, the diameter of rotor is 120 
meters and 100 meters height for the tower. 

Wind farms are generally located on the area of strong 
wind speed such as on mountain areas. They are 
usually called land-based wind turbine farms. Due to 
the higher expected energy, wind farms are moved 


from land to offshore as the wind speed is higher there. 
In order to anlyze wind turbine in offshore condition, 
the coupling of different physical fields are necessary. 
Here we present how to model such kind of the 
problem. 

Entire Model of Offshore Wind Turbine 

The model of an offshore wind turbine consists of four 
sub-models, the aerodynamics of wind, the structure 
dynamics of turbine components, the hydrodynamics 
of the ocean wave and the soil dynamics of 
unbounded soil, see in Fig. (1). All sub-models are 
developed initially separately and finally coupled via 
the partitioned method. The descriptions of each sub- 
model are as follow; the structural subsystem, which is 
a nonlinear large displacement finite element model 
(FEM); the aerodynamic subsystem which describes 
stochastic wind and includes instationary dynamic 
stall model for the aeroelastic blade loading; the 
hydrodynamic subsystem, modelled mainly as 
potential flow and couple to the stochastic wave field 
to describe the real sea state, the computations on this 
part based on the boundary element method (BEM); 
and the soil dynamic subsystem modelled as near- 
field soil near the underground structure discretized 
by FEM and far-field effects included by the scaled 
boundary finite element method (SBFEM). 


wind 



FIG. 1 ENTIRE MODEL OF OFFSHORE WIND TURBINE 

All of these models are coupled in time domain, and 
the consistent simulation of the offshore wind turbine 
can be performed based on partitioned method. 
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Structural Dynamics of Turbine 

On the structural dynamics part, we followed 
previous research from our institute. The geometrical 
nonlinear formulations of the beam are used for 
dynamics of structure. The structure of turbine is 
initially seperately modeled in substructure level, 
three blades structure for rotor system, tower system 
and foundation system. 

Formulation of Coupled Substructure System 


two steps; the first is to transform from the rotor 
system to the intermediate reference system (index a) 
at the coupling node, while the second is to transform 
from the intermediate reference to the tower system 
(system i = 0), N ° at the coupling node, the 

transformation can be defined by the following matrix 
~N% 0 O' 

N% = 0 W?(2,2) W?(2,3) , (5) 

.0 JV?(3 # 2) JV?(3 # 3). 

with coupled variable qf = [u 1 ,u 2 ,u 3 ,^ 2 ^ 3 ] r • 


First of all, we focus on the rotor system. The rotor has 
three blades. Each blade has the variable q l associated 
with it, where (i = 1, 2, 3) denoting the three blades. 
The variable q l is represented by a coupling variable 
(subscript c) and the free variables (subscript f) as 
shown below 



Similarly, the tangential stiffness D l and the force 
vector g l for each blade are as follow 



(2) 


g l 


9c 

.aY 


(3) 


The blades are rigidly coupled to each other, each 
blade has its own coordinate system. We choose one 
blade (i = 1) as the reference coordinate system for the 
rotor system. The transformation of this reference 
system into other two blades (i=2,3) can be performed 
with the transformation matrix N\ and N\ . The 
variables are given by q 1 = N\q 2 and q 1 = N\q 3 , 
where N is transformation matrix, see more details in 
reference. After rewrite the formulation above yields 
the coupled system of rotor blade. 

D r q r = (4) 



Now the rotor system has to be coupled to the tower 
system (numbered i = 0) according to the assumption 
that the rotor blade can freely rotate about the axis (the 
axis of wind direction). The axis of the wind direction 
may be defined as the X-direction of the coordinate 
system (i = 0) of the tower. During the computation, 
the transformation from the rotor coordinate system 
(system i = 1) to the tower coordinate system (system i 
= 0) is required. This transformation can be done in 


Finally the coupled system of structures can be 
achieved by substituting and coupling the variables as 
shown below 

Dq = g, (6) 

D ff D % 0 1 r q0f 

where D= D° cf D° cc + N\D r cc N? N\D r cf ,q= q° c 
_ 0 D} c Nf D r ff J W/ 

g° f 

9= g°c+Nl9 r c • 

i 9 } J 

The system of equation above shows the coupled 
substructure system of turbine components. The same 
procedure is performed to complete the dynamics 
equation, in which the effect of rigid body rotation of 
rotor system is taken into account. 

Dynamics Equation ofTurbined 

The dynamics equation of wind turbine is described 
by system of partial diferential equation (PDE). The 
formulations are derived based on geometrical 
nonlinear of beam, the systems of PDE are shown very 
brieftly here, more details please refer to reference. 
The dynamics equation is written in two balance of the 
momentum, forces due to translation and moment due 
to rotation. 

The balance of forces is 

d,(J?n e ) + n ext = ndfu, (7) 

where R is rotational matrix, n is forces. The balance 
of moment given by 

d^Rm e ) + d,r x (Rn e ) + m ext = R{j]Y ) +Jd t Y,(8) 

The Y is the skew symmetric matrix of angular 
velocity Y, and J is inertia tensor. The Eq.(7) and Eq.(8) 
are discretized based on conventional FEM. The 
discrete form of dynamics equation is derived and 
yields as below 

f k (.dtd, d t d, d, t) + fi n (d, t) = f ext (d, t), (9) 

where f k denotes inertia forces, in each element the 
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inertial force given by 

fe _ C l k+ 1 

Tk ) lk 


l-idtd 

T T (Jd t Y + d L Y]Y. 


the internal forces in element can be computed 

\ C nY] 


r l k+l [ T r 

J ik L N J L< 


C m K. 

where v denotes corresponding matrix 
[C n y C m K] T denotes corresponding matrix of material 
The external forces f ext is interaction force from 
coupled external system, e.g. wind, wave, sea bed. 
This forces given by 



Wind 



(10) 

fe 



T 

JL Wake 

(11) 



Structure 

and 





f, 


l k + 1 

Ik 


N 1 


RHext 

T T R T m p 


(12) 


where n ext denotes external forces and m ext denotes 
external moment. 


I 


FIG. 3 THE CONFIGURATION OF FLOW FOR WIND TURBINE 
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Aerodynamics of Instationary Wind 


Rotor position 


On dynamics of wind turbine, the aerodynamics of 
instationary wind are required for computing wind 
forces along the rotor blade. This part we followed 
previous research of our institute. The model 
described local fluid flow through profile of rotor 
blade. The global flow is also necessary in order to 
include the wake effect on rotor system. 

LEft 

Rotations 
rtiavemetiL 

Profik- 
-vihrjdkin 


Wind 

FIG. 2 FORCES ON 2D PROFILE OF ROTOR BLADE 
Local Aerodynamics 

Fluid flow model describe the flow through rotor 
blade locally. The forces along 2D profile of rotor 
blade are reasonably obtained in this case, see in Fig. 
(2). The drag force (f d ) acting parallel to the relative 
velocity ( v r ), while lift force (f t ) is perpendicular, the 
moment ( [m ae ) applies on 2D profile for take the 
rotational effect due to local aerodynamics flow. The 
formulations are given by 


fd =e f v r 2 C d C, 

(13) 

fl =e f v r Z Cl C ’ 

(14) 


(15) 


where p a denotes air density, c is cord length of profile, 
while c d ,c t , c m are drag, lift and moment coefficient, 
respectively. The relative velocity ( v r ) is the means of 
wind velocity, this ( v r ) includes the effect of rotational 
movement of rotor and the vibration of rotor blade. 




FIG. 4 THE IDEAL FLOW OF ROTOR SYSTEM 

Global Aerodynamics 

As the global aerodymics of wind turbine is very 
complicated due to turbulence and vortex behind the 
rotor (see in Fig. (3)). These physical flow could be 
modelled and solved them somehow by CFD, 
however the computational effort is very limited. 
Here we used the ideal flow of rotor system, under the 
assumption that no flow in radial direction inside 
wake and flow is frictionless. The wind speed is 
reduced from upstream with velocity v 0 to down 
stream with velocity v 1 . The pressure increase from 
atmospheric pressure at upsteam position to pressure 
p at rotor position before it is immediately dropped A p 
by the rotor, see in Fig. (4). Again at downstream, the 
pressure increases to the atmosmeric level. The acting 
forces on rotor / rot can be obtained by 

frot ~ ^pA r ot> (16) 

where A rot denotes area of rotor dish. 

The pressure drop A p can be computed by using 
Bernoulli's equation from upstream to downstream, 
which yields the pressure dropp as 

Ap = - 2 p{y o 2 (17) 



FIG. 5 THE CONTROL VOLUME OF WIND FLOW 
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Concerning the control volume in Fig. (5), the axial 
momentum applied to this volume, yields the 
equation below 

pv 1 2 A 1 + pv 0 2 (A c - AJ + d t m.v 0 - pv 0 2 A c = -F,( 18) 

where d t m is substituted by means of mass 
conservation. The forces F coming from pressure drop 
over the rotor, when we substitute force (F) in Eq.(18) 
by pressure drop from Eq.(17), yield the interesting 
formulation. 

V = ^(v 0 + v 1 ). (19) 

Eq.(19) shows the physical meaning that the velocity at 
rotor position can be obtained by means of velocity at 
upstream and downstream. This makes less 
computation to obtain the velocity at the rotor blade. 

Until now the computing aerodynamics force to rotor 
blade has presented. These are so-called blade element 
momentum theory, however, the formulations are 
suitable for infinite number of rotor blade. When used 
for finite number of blade, the correcting factor is 
necessary, which can be obtained by the Plandlt's tip 
loss factor, more details please refer to reference. Also 
the instationary of wind can be involved by using 
dynamics stall model. The systems of ODE are 
presented to describe the instationary effect of 
aerodynamics forces as shown in Fig. ( 6 ), and there is 
unnecessary here to go in details, readers please refer 
to [3] and reference therein. 
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FIG. 6 INSTATIONARY O F AERODYNAMICS COEFFICIENT 


In order to complete aerodynamics forces on the blade 
structure corresponding to FEM in structure part, the 
aerodynamics forces along the element of the rotor are 
computed as 


fe _ C l k+ l 
Ta ~ J lk 


n t 


R^ext 

T T R T m e J 


(20) 


where R and T are corresponding matrix in FEM. The 
external forces n ext = [ 0 , cy/ t , Cff n ] and external 
moment m ext = [cym ae , 0 , 0 ] , where cy denotes 
PrandtTs tip loss factor. To complete f e a in Eq.(20), 
there is suitable to integrate numerically. The Grauss 
point integration scheme is used here, as this reqires 
only the informations of n ext and m ext at Grauss point 


in 2D case, see in Fig. (7). The coupled aerodynamics 
of wind and dynamics of structure is automatically 
included in the model as mentioned before, in 
addition, the aerodynamics forces ( n ext and m ext ) are 
function of relative velocity v r and aerodymics 
coefficient (c d ,Ci,c m ). Both are function of displace- 
ment , velocity, and acceleration of structure. These 
mean that when the movement of structure is changed, 
it influence to wind flow and vortex system. 
Oppositely when the wind flow is changed, it 
influences the movement of structure as well. 
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FIG. 7 AERODYNAMICS FORCES AT GRAUSS POINT 

Hydrodynamics of Ocean Wave 

On the part of ocean wave, we model the fluid flow 
based on potential theory. The fluid is incompressible, 
inviscid, irrotational flow and the surface tension is 
nelegted. The velocity field is obtained from gradient 
of potential O as v(r, t) = VO in any position (r) and 
time (t) in fluid domain (Qf), where V denotes the 
spatial gradient reference to the corresponding 
coordinate system. The fluid inside domain 12 f has to 
be satisfied by so-called Lapace's equation. 

A 0>(r, t) := V 2 0(r, t) = 0. (21) 

The free surface of fluid is described as level set 
function F(r, t). The fluid particles on the free surface 
are moved together with the free surface movement. 
The normal velocity of fluid can be computed by 

vn = V4> n =: 3 n V 2 <l>(r, t) = (22) 

where n is normal velocity of fluid on surface. 

On the free surface, it will be convenient to introduce 
the coordinate system, in which z-direction is opposite 
to gravity and point upward from undisturbed x,y 
plane. The level set function is obtained by following 
formulation 

F(r, t) = F{x , y , z, t) = r](x, y, t) - z = 0. (23) 

The kinematic free surfce boundary condition (KFSBC) 
shows no flow through the free surface. 

D t F(r, t) = d t F - VO ■ VF = 0, (24) 

concerning the level function on x,y plane (at z= 0 ), the 
free surface elevation 77 satisfies KFSBC 
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d t r] = d z <& - V h 0> ■ V h rj = 0, 25) 

where V h denotes spatial gradient of horizontal x,y 
plane. The dynamic free surface boundary condition 
(DFSBC) is given by Bernoulli's equation. 

d t <P =gt] + 1 -\V<t>\ 2 -^ = 0, (26) 

z Pw 

where g is gravity, p a denotes pressure and p w denotes 
density of fluid. Similar to the case of free surface, the 
kinematic and dynamics boundary condition on the 
fluid structure interface can be written as in Eq. (25,26). 
Only the particular coordinate system should be 
introduced, the g-axis is normal to the plane (4,$) . 
The movement of the fluid structure interface is 
introduced as x(t 0 = Xo g, t) , 

where x 0 (£,,$,g) denotes the position of interface at 
rest and g,t) denotes the displacement of the 

interface. Let Xqi^> Q 0 is the movement of interface 
in g -direction. The kinematic boundary condition 
describes no flow through the interface 

d t X q = - V,0 ■ V s * f , (27) 

where V s denotes surface gradient on (£,$) coordinate. 
On the boundary condition of sea bed, the d t x<; = 0 
describes no movement of sea bed, oppositely, if the 
movement of sea bed is concerned, the Eq.(27) can be 
used for the interface such as zunami wave. The 
dynamics boundary condition is given by Bernoulli's 
equation. Only the interaction forces between fluid 
and structure have to be included, containing stress 
vector . 

=g(z + r 1 + u z ) + 1 - |VO| 2 -^,(28) 

z Pw 

where u z denotes z-component of displacement in 
global coordinate system (x,y,z). 

The random finte depth wave is also developed in this 
research work, in order to make the possibility of 
including random wave and fully nonlinear wave into 
the computation. The details of mathematical model 
will not be shown here again, however these have 
already been published and are available in references. 

Discretized Formulations with BEM 

The Laplace's equation is initially written in the 
integral form by using Green's second identity, which 
yields the boundary integral equation [6,7]. 

C(x r ) 0(* r ) + J &(x) d n G(x,x r )dr(x) 

- f d n O(x ) G(x,x r )dT(x ) = 0, (29) 

where x r (x r ,y r ,z r ) denotes the refrence point and 
x(x, y, z) denotes the arbitary point on the boundary. 
The term G (, x r ) is obtained by 


C(* r ) = £ ( 30 ) 

where a is the solid angle on the boundary, to 
compute solid angle, please refer to reference. The 
corresponding Green's function for 3D case is given by 

l 


G(x, x r ) 


4n\x-x r \’ 


(31) 


The BIE in Eq.(29) is discretized by BEM, and meshes 
are necessary only on the surface (2D) of boundary for 
3D domain of fluid. The variables on the boundary 
node are O and d n O , in each node one variable should 
be defined to complete the boundary value problem. 

A O = Bd n O, (32) 


Eq.(32) above is rewriten by setting known value p k 
on the right hand side and unknow value p u on the 
left hand side with corresponding matrix H and G . 

Hp u = Gp k . (33) 


In order to reduce the CPU time, the Fast Multipole 
method is also implemented here, more details about 
Fast Multipole method for BEM please refer to 
reference and reference therein. 


Soil Dynamics of Sea Bed 

Mathematical Formulations of Sea Bed 

The formulations of soil dynamics for sea bed will be 
explained here brieftly, for more details please refer to 
[6,9]. The sea bed is modelled as continua based on 
linear elastic theory. The stress (a) and strain (e) field 
are represented by Hooke's law. The constitutive 
relation can be written as 

g = Ec, (34) 

where E is elasticity matrix, the strain field can be 
obtained by strain-displacement relation. 

e = Pd, (35) 

The differential operator P is combined in cartesian 
coordinate system. Based on Newton's law, the 
dynamics equilibrium equation in frequency domain 
can be obtained as 

P T o + b = -a) 2 d, (36) 

while Eq.(36) can also be written in time domain as 

P T o{f) + b{t) = pd}dit), (37) 

where stress a, body forces b and displacement d are 
in time domain. 

Discretized Formulations with FEM 

The dynamics equilibrium equation is dicretized based 
on the conventional finite element methods. Here the 
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sea bed is considered into two regions, in which a 
finite region represents soil close to the structure, 
including some parts of structure (near field) and an 
infinite region represents unbounded soil far from the 
structure (far field). The discrete form of Eq.(37) is 
given by 

M c d?d(t) + K c d(t) = f c (t), (38) 

where M denotes mass matrix and K denotes stiffness 
matrix while / is force vector given by 


-Mss 

M s b 

,K C = 

Kss K sb 

-Mbs 

M bb - yMM°Z. 

_K bs Kbb. 


f _r fs(p: 

) 1 



The subscript c denotes the coupled system of near 
field/far field. The subscript s and b are inside domain 
and at the boundary respectively. The effect of infinite 
domain will be included in yAtM™ term with y 
paramter in Hilber- Hughes- Taylor (HHT) time 
integration scheme. The r h (t) denotes interaction 
forces from infinite domain. 



\ infinite domain 

FIG. 8 FINITE AND INFINITE OF SEA BED 


Scaling Center Near field / far field interface 

FIG. 9 SCALED BOUNDARY TRANFORMATION 

Coupled with SBFEM for Unbounded Soil 

In order to include the effect of unbounded soil, the 
coupled near field/far field is developed as mentioned 
in previous section. The far field is involved here by 
using the so-called Scaled Boundary Finite Element 
Method (SBFEM), and the concept of geometrical 
similarity is implemented associated with conventional 
approaches of FEM. Initially the near field soil 
(including some part of srtucture) must be discretized, 
FE-meshes will be constructed in this part, while 
SBFE-meshes can be obtained by using FE-meshes on 
the boundary, see in Fig. (8). The FE-meshes (on the 
boundary) are transformed into the so-called scaled 
boundray coordinate system (rj, (, f), see in Fig. (9). 
The axis rj, ( is described in the circumferential 


direction only on the boundary, while the third axis ( 
contains a scaling factor from center. This scaled 
boundary coordinate permits a numerical treatment in 
circumferential direction rj, ( based on weighted 
residual technique. The PDE in Eq.(36) is transformed 
into ODE in radial coordinate ,yield as 

E 0 (f 2 d^d + 2 - E^d^d + E\{^d^d + d) 
-E 2 d + co 2 M 0 f 2 d + ff t +efb = 0, (39) 

where E 0 ,E lf E 2 ,M 0 are the coefficient matrix which 
are independent of f. The f t ,f b denote traction forces 
and body forces. 

Continously, the method of vitual work is used to 
derive the SBFEM formulation in term of dynamics 
stiffness 5°°of infinite domain. 

(f^S 00 + EJE^^S 00 + El) - f^S 00 - 
fdf(C 1 S~)-E 2 +Z 2 M 0 = 0. (40) 

The Eq.(40) still in frequency domain (5°° = S°°(co)). In 
order to transform Eq.(40) into time domain, the 
inverse Fourier transformation is used to derive the 
formulation in time domain. It is not necessary to go 
in details of deriving of formulation here, autors refer 
to references. Finally the acceleration of unit impulse 
response matrix M°°(t) can be obtained by the 
following formulation 

M°°(t) = U T M°°(t)U f (41) 

where M°°(t) can be determined by solution of 
convolution integral of dynamics formulation. U T U are 
decomposition of coefficient matrix E 0 . 

In Eq.(38), the interaction forces r b (t) at near field/far 
field interface are required which could be obtained by 
means of convolution integral, here we show the 
discrete form of these formulation 

n 

r h (t) = ^ M~_j J dfd(T) dT 
J= 1 0'-i)a t 

= Yj=i MZ-j (d t dj - dfdj-j). (42) 

For large number of n (i.e., long time simulation) the 
direct solution might be limited due to the very time 
consuming. The recursive procedure is used for 
approximation this forces with reduced computing 
time, in details of recursive procedure the autor refer 
to reference and reference therein. 

Coupled Multi-physics Model 

The general pattern is that in each subsystem the 
dynamics are described by either pure differential 
equations (abstract ODE), such as for the structure 
Eq.(9), or by differential algebraic equations (abstract 



6 


Development and Applications of Oceanic Engineering (DAOE) Volume 3, 2014 


www.daoe-journal.org 


DAE). The total coupled system as described before by 
the space-discretized version may be written as the 
following form 

d t <P = (43) 

dtV = <&,(<*>, */,T,d), (44) 

0 = ’/'<!>, i; (*, V’ T > d), (45) 

M(d)dtd + ft nt (.d t d, d) = 4> s (d t d, d, <t>,?j),(46) 
0 = ^ (<P,V.r,d), (47) 

M{d)dfd + fint(d t d, d) = <j> h (d t d, d, <J>,j/),( 48) 

0 = 'F ssi (r,d), (49) 

d t x = 4> w (d,d t d), (50) 

0 = d, d t d). (51) 

The Eq.(43)-Eq.(51) describe all subsystem of modeling, 
however, the modeling of subsystems can be changed 
without difficulty, i.e. fluid modeling can be changed 
from incompressible to compressible, while the 
coupling procedure stays the same but only subsystem 
is changed. To so vie them in time domain, they have 
to be solved based on iterative coupling procedure 
(partitioned method). The coupled algorithms are 
described here with the general flow of information. 
The subsystems are separated softwares, so called 
software components, coupled via the Component 
Template Library (CTL) used as the middle- ware in 
order to communicate with all the different software 
components. The coupled algorithms are described in 
the following simple way, see also Fig. (10).: 


A similar coupling exists between structure and soil. 
The structure passes displacements (Dirichlet data) to 
the soil through the soil-structure interface. The 
dynamics of soil are solved with the HHT method and 
forces (Neumann data) are returned to the structure. 
The structure gets the corresponding forces from the 
coupling interface from soil and fluid. A new 
Newmark trial-step can be computed. 

Also similarly coupled with the aerodynamics and 
structure, the dirichlet data is given from the structure 
part at the coupling area (here it is on the blade). The 
aerodynamics are solved and compute the 
aerodynamics forces. It returns the aerodynamic forces 
as the Nuemann data to the structure component to 
complete the structure dynamics. The iteration loops 
can continue in parallel until the residuum is 
sufficiently small. 

Software Component 

As in the coupling multi-physic model, there are many 
physical domains that have to be coupled. In each 
model, the corresponding software has to be used as 
the component of coupled systems. They are 
implemented in the different ways. The structural 
dynamics subsystem, the nonlinear large displacement 
modelling is used and discretized with FEM. The main 
code of this component is implemented in MatLab and 
some parts which take more CPU time are embeded 
with C functions. It works as the normal command 
line in MatLab, but it runs on C code. 


« 

>>/ 5v/ 

S IR! r< 

■y/, "x\ 



FIG. 10 COULING PROCEDURE FOR ENTIRE SYSTEM 

At the current time step, the hydrodynamics quantities 
are computed with the PC-algorithm, termed as 
predicted values. This involves the solutions of 
Laplace's equation for each of the fluid domains. Now 
the hydrodynamic forces can be computed at the fluid- 
structure interface, and they are passed to the 
structure at the current time (Neumann data ( labeled 
N). The structure system is integrated with these 
forces. The computed displacements change the fluid 
domain, the structural velocities are Dirichlet data for 
the fluid domains (labelled D) in Fig. (10) when 
passing back to the fluid. The iteration then goes back 
to the first step. This describes the fluid-structure 
coupling. 


The aerodynamics subsystem, the blade element 
momentum theory is used combine the stochatic wind 
field. This part is implemented in MatLab and 
communication with the structural part with the 
Matlab function. The hydrodynamics subsystem as the 
potential flow with the stochastic wave field, the main 
code of this is implemented in MatLab. Just only the 
part of solving Laplace's equation, we implemented 
the fast multipole-code FastLap (it is written in 
Fortran). To communicate Matlab code, the embedded 
C-f unctions is also used here. 

The soil dynamic subsystem is modelled as the near- 
field and far-field model. The near-field is discretized 
by FEM, while the far-field is discretized with the so- 
called scaled boundary finite element method 
(SBFEM). The FEM is implemented with the Felt code 
(it is written in C) and the SBFEM is implemented with 
the similar code (it is written in Fortran) . 

For the realization of this coupled simulation, 
graphically depicted in Fig. (11), we have chosen the 
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concept of software-components. Here a component is 
a set of software which consists of a well defined 
interface and an implementation. Interface and imple- 
mentation are connected through a communication 
channel. In this way, the usage of a component is 
independent of its location. 



FIG. 11 SOFTWARE COMPONENT 


The middle-ware used to implement the component 
technology is the Component Template Library (CTL) 
based on C++ generic template programming. It can be 
used to realise distributed component-based software 
systems. This library serves as an easy-to-use 
programming environment for distributed applica- 
tions in an abstract manner, but its main focus is to 
transform existing C/C++ or FORTRAN libraries to 
remotely accessible software components. As underlying 
communication protocols MPI, PVM, or directly 
sockets, as well as dynamic linkage and threading are 
supported. 

Simulation Results 

Firstly, we show the numerical testing of soil 
subsystem coupling with tower (no effect rotor) in 
order to study the unbounded soil model with 
conpled FEM/SBFEM as explained in section of soil 
dynamics of sea bed. 


Applied forces 



The monopile is located at the depth d2 = 20 m into 
finite domain of the soil, in which the depth of the soil 
domain h = 30 m. and domain width w = 20 m.. The 
length of the tower is measured from the surface of the 
soil dl = 80 m. Horzontal force is applied on the top of 
tower. The dynamics equation of the soil is solved 


with the HHT method combine with the constant time 
interval 0.005 sec. 

The result shows effect of unbounded soil on the 
movement of structure, see Fig. (13). The amplitude of 
displacement decreases along the time, meaning that 
the energy is damped due to unbounded soil, which is 
very physically reasonable. Therefore this FEM/ 
SBFEM is suitable for our entire simulation of off shore 
wind turbine. For other couple subsystem such as 
wind-structure or wave structure, we have also tested 
the coupled subsystem separately and these have been 
published. The entire model of offshore wind turbine 
have simulated based on the dynamics coupling 
procedure as described before. Concernig the 
simulation result, on the part of hydrodynamics, the 
wave formulations are written based on the fully 
nonlinear assumption as described before. The 
random wave fields are included by the stochastic 
process at the far from the structure. These waves are 
actually linear and the irregular wave charateristics 
are allowed to be involved as the summation of many 
differrent waves. The waves propagate and develop 
itself with the fully nonlinear formulation. The waves 
became nonlinear and its steepness increased during 
the propagation. 


E. 

S 

1 
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FIG. 13 TIME HISTORY OF DISPLACEMENT AT TOP OF TOWER 

The irregular waves form remains during this 
propagation. When wave arrived the position close to 
structure, these waves were developed to become 
nonlinear and irregular waves. The waves run through 
the structure and pass away, while some waves are 
reflected by the effect of structure vibration. 


o 

! 
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FIG. 14 IRREGULAR WAVE NEAR THE STRUCTURE 



Number of time step 
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FIG. 15 VISUALIZATION OF WAVE INTEACTING STRUCTURE 

Fig. (15) shown the nonlinear wave interacting the 
tower, and the reflection wave can be observed due to 
structure vibration. Also the irregular wave form can 
also be found consequently as shown in Fig. (14). 

Considering the result on rotor blade and tower, the 
result shows the irregular form of deformation, 
resulting from irregular form of wind and wave. The 
deformation of rotor blade and tower are shown in Fig. 
(17), while the displacement of rotor is shown in Fig. 
(16). 



FIG. 17 VISUALIZATION OF BLADE AND TOWER 



FIG. 18 VISUALIZATION OF OFFSHORE WIND TURBINE 


Fig. (18) shows the visualization of entire system of 
offshore wind turbine. The deformation rotor blade 
and tower can be observed as their movements are 
very irregular form. In front of rotor, we present the 
velocity of the wind field, and this wind velocity is 
very random and presents irregular form of wind 
forces acting rotor blade. 

Conclusions 

In this paper, the model of an offshore wind turbine 
was presented consisting of the aerodynamics of the 
wind, structure dynamics of the turbine component, 
the hydrodynamics of the ocean wave and the soil 
dynamics of the sea bed. A concept of modelling 
offshore wind turbines is to develop the sub-system 
modelling initially separately, and then coupling the 
entire model by using partitioned method. The 
structure modelling of offshore wind turbines is 
modelled by the nonlinear beam element and 
discretized with FEM. The random of wind field is 
included based on the blade element momentum 
theory and dynamics stall model. The hydrodynamics 
of the ocean wave is modelled as the random finite 
depth wave. The fluid dynamics is described by the 
BIE of Laplace's equation and discretised by BEM. 
The soil can be modelled by the finite domain and 
infinite domain. The finite domain is discretized by 
finite FEM, and the infinite domain can be included by 
SBFEM. 

The software components are communicated with 
CTL middle ware. The simulation of coupled multi- 
physics model was performed in time domain and 
reasonable result was presented physically. However 
in either sub-modelling the complexity could also be 
extended or increased in order to come closer to the 
nature, which would be outlook of our research and 
will be presentedlater. 
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